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General relativistic polarized radiative transfer: 
building a dynamics-observations interface 
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ABSTRACT 

The rising amount of polarized observations of relativistic sources requires the correct 
theory for proper model fitting. The equations for general relativistic (GR) polarized 
radiative transfer are derived starting from the Boltzmann equation and basic ideas of 
general relativity. The derivation is aimed at providing a practical guide to reproducing 
the synchrotron part of radio & sub-mm emission from low luminosity active galactic 
nuclei (LLAGNs), in particular Sgr A*, and jets. The recipe for fast exact calculation 
of cyclo-synchrotron emissivities, absorptivities, Faraday rotation and conversion co- 
efficients is given for isotropic particle distributions. The multitude of physical effects 
influencing simulated spectrum is discussed. The application of the prescribed tech- 
nique is necessary to determine the black hole (BH) spin in LLAGNs, constraining 
it with all observations of total flux, linear and circular polarization fractions, and 
electric vector position angle as functions of the observed frequency. 

Key words: black hole physics - Galaxy: centre - plasmas - polarization - radiation 
mechanisms: general - radiative transfer 



1 INTRODUCTION 

The good model of radiative transfer is the key in bridg- 
ing the plasma dynamics and the observations of com- 
pact accreting sources . The dynami cs of plasma evolved 
from hydrodynamics {Ruff ertl Il994f) to m agneto hydro- 



dynamics (MHD) (lHawlev fe Balbusi 120021) and p article - 
in-cell (PIC) simulations (|Sironi fe Spitkovskvl |2009| ). 
The modelling of compact objec t's gravity has t urned 
from quasi-Newtonian potentia l (JHawlev fe Krolikl l200ll ; 
llgumenshchev fc Narava n 2002) to the full general rel- 
ativistic (GR) MH D (|De Villiers. Hawlev fc Krolikl 120031 ; 
IShafee et all |200S| ). Only GRMHD simulations allow to 
fully account for the spin of the compact object. The ra- 
diative transfer approximations were improvi ng as well. 
The simple quasi-Newtonian ray propagation (|Chan et al.l 
120091) gave way to null-geodesies tracing in Kerr met- 
ric (ISchnittman. Krolik fc Hawlevl l200rj ; iNoble et all 120071 : 
iMoscibrodzka et al.ll2009i ). The huge amount of polarization 
observations demanded the polarized radiative transfer. 

The main principles of GR po larized radiative transfer 
were formulated in 1 Broderickj (2004). However, that formula- 
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tion was not ready for applications as it lacked, for example, 
the Faraday conversion and the suppression of Faraday rota- 
tion in hot plasmas. The first a pplication to the real object 
was done in lHuang et alj (|2009i ). Their calculations included 
Faraday conversion, but made several approximations, some 
of which can be substantially improved upon. For example, 
their simple relation on V and Q emissivities and Faraday 
rotation and conversion constitutes an approximation that 
almost never holds. Their emissivities are calculated in syn- 
chrotron regime, which breaks for temperatures about the 
electron mass. Their frame of plasma does not fully account 
for the fluid motion. We are improving on their work in 
the present paper, in particular treating exactly the plasma 
response and extending it to non-thermal particle distribu- 
tions. 

Another important issue is the complexity of GR po- 
larized radiative transfer. The errors and implicit strong ap- 
proximations may slip into the equations of almost every 
author. This is more likely the case, when certain deriva- 
tion is done half-way by one author and then continued by 
another author, e.g. the derivation of the Faraday conver- 
sion coeffici ent in the mix ture o f thermal and non-thermal 
plasm as in I Melrose! i|1997f ) and iBallantvne. Ozel fc Psaltisl 
(|2007l ) neglected the importance of the finite ratio of cy- 
clotron to observed frequencies. Another good example is 
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right-handed rotation along the ray: 6>0 



Figure 1. Geometry of the problem. Vector Bo represents uni- 
form magnetic field. The transverse plane wave travels along fc 
and has electric field E oscillating in (e^-e 2 ) plane. Vectors a and 
b represent parallel transported basis orthogonalized with fc. 



the definition of the sign of circular polarization V. It varies 
from article to article and the consistent definition in a single 
derivation is essential. 

Therefore, there appeared a need for the present paper. 
In a single derivation from the basic principles we provide 
the necessary applied expressions for GR polarized radia- 
tive transfer. We start in 32]by consistently defining the po- 
larization tensor, Stokes parameters, and plasma response 
tensor and incorporating the response tensor into the New- 
tonian radiative transfer. In [f3] we recast the derivation of 
the response tensor from Boltzmann equation for general 
isotropic particle distribution and do the special case of 
thermal distribution. We provide the applied expressions 
for response tensor in the plane perpendicular to the ray 
and give the consistent sign notation for both positive and 
negative charges in H3.3I The resultant formulas for absorp- 
tivities/emissivities, Faraday rotation and conversion coeffi- 
cients are exact and easy to evaluate. By means of locally- 
flat co-moving reference frame we extend the radiative trans- 
fer to full GR in 33] We highlight the various physical effects 
important for real astrophysical objects in Sj5] Finally, in [J6] 
we briefly summarize the methods and the ways to general- 
ize them even further. 



2 NEWTONIAN POLARIZED RADIATIVE 
TRANSFER 

The proper treatment of polarization of radiation is neces- 
sary to take the full advantage of polarized observations. Let 
us start formulating the dynamics of polarization by defin- 
ing the basis. Let e 3 be the direction of uniform magnetic 
field Bo. Then def ine the orthonormal triad fc, e 1 , and e 2 in 
the standard way (jRvbicki fc Lightmanll 1979 ; ISazonovll 1969 : 




Figure 2. Right-handed rotation of electric field along the ray at 
fixed time t corresponding to negative circular polarized wave V < 
0. Electric field is E\ = E\^ exp(i(kz + 8)), E 2 = E 2u} exp(ikz), 
8 = 7r/2, where Bi^ is along e 1 , E 2uJ along e 2 , and the wave 
propagates along fc. Vectors e 1 , e 2 , fc constitute the right-handed 
triad. 



IPacholczvkll 197(3 ): wave propagates along fc vector, 

e 1 = C(B x fc), 

2 , 1 

e = fc x e , 



(1) 



where the scalar C can have either sign. We choose the axes 
as on Fig. [T] e 1 is perpendicular to (Bo,k) plane and Bo 
lies in (fc, e 2 ) plane. The rotation around e 1 transforms basis 

e 3 ) as 



(e\ e 2 , e 3 




(2) 



(3) 



Vectors and tensors then rotate according to 

A k = (M T ) k iA t , a ki = (M T ) fcm Q m „M„,, (4) 

where () is a transposed matrix and the quantities with 
tildes are taken in a frame with the basis (e 1 , e 2 , e 3 ). The 
angle 6 can be found from 

„ fc • Bo 



e 3 = fc = e sin 9 + e 3 cos 9, 
which can be conveniently written as 

e fc = e l M lfe , 



kB 



(5) 



For the electric field components E\ along e 1 and E2 
along e 2 the Stokes parameters are defined as 

/ = {EiEt) + (E 2 E* 2 ) , 

Q = (EiEj) - (E2EI) , (6) 

U = {EiE£) + {E 2 E{) , 
V = i((E 1 EZ)-(E 2 El}), 
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clockwise rotation in fixed plane: <5>0 




Figure 3. Left-handed rotation of electric field at fixed co- 
ordinate z with time corresponding to negative circular polar- 
ized wave V < 0. Electric field is E\ = E\^ exp(«(— uit + 8)), 
E 2 = E 2lj exp(—iu}t), S = tt/2, where Ei^ is along e , E 2ul along 



where the last formula cho oses the IAU/IEEE definition 
IjHamaker fc Bregmanl 1 19961 ) of V, actively used by ob- 
servers. That is for positive V the rotation of electric field 
is counter-clockwise as seen by the observer. Nevertheless, 
all ast rophysics textbooks agree o n the opposite definition 
of V dSazonov fc TsvtovIcbi1l968l; iLegg fc Westfoldl ll96Sl; 



Rvbicki fc Lightmanl Il979l : iRochfordl T2001I : IWilson et all 
20091 '). Let us visualize the electric field rotation. Take the 
electric field 

Ex — Ei^exp^kz-ujt + S)), E 2 — E 2u; exp(t(kz-oJt))(7) 

with positive amplitudes (Fourier coefficients) E\^ , E 2uj > 
and substitute it to the definition ©. Then 



+ eL 



I = E( u 

U = 2E luJ E 2ul cosS, 



Q = E lhJ — E 2uj 
V = -2E luj E 2uj smS 



(8) 



and V < for S e (0,7r). Let us fix time t — 0, 6 — ty/2 
and draw the electric field vector in space along the ray (see 
Fig- E| . We see that the electric field corresponds to a right- 
handed screw. However, if we fix a plane in space by setting 
z — and draw the evolution of the electric field vector, then 
the rotation direction is the opposite: electric field vector ro- 
tates counter-clockwise, if viewed along the ray (see Fig. [3]). 
These opp osite directions of rotation is a common point of 
confusion (|Rochfordl l200lh . Correspondingly, the observer 
sees the clockwise-rotating electric field for 5 — n/2, as she 
is situated at a fixed z, and counter-clockwise rotating elec- 
tric field for 8 = —tt/2 and positive V > 0. The definition 
$6$ of V has a marginal advantage: the electrons generate 
positive V for propagation along the magnetic field. Let us 
take an electron on a circular orbit in (e 1 ^ 2 ) plane. It moves 
from the direction of e 1 to the direction of e 2 , id est clock- 
wise as viewed along the magnetic field. Then the wake of 



the electric field follows the charge and rotates clockwise. 
The resultant wave propagates along Bo and constitutes a 
left-handed screw, which gives the positive V > 0. 

The polarization tensor is obtained automatically from 
equation (0 as 



hi = (EiE*) = I ( 



I + Q U-iV 
U + iV I-Q 



(9) 



The polarization vector is 

S=(I,Q,U,V) T . (10) 

Note that iMelrose fc McPhedranl (|l991f ) uses the same def- 
initions, however, their fc, e , e 2 constitute a left triad in- 
stead of a right triad. 

The plasma response is characterized by the 4 x 4 re- 
sponse tensor a M „ 



./: 



o? v K 



(ii) 



as the proportionality between the four- vectors of vector po- 
tential amplitude and the current density amplitude. There 
is freedom in choosing the gauge condition for A M . Let us 
choose the Lorenz gauge 



A»u u 







(12) 



at each point along the ray and enforce it by adding to A^ 
a vector, proportion al to fc^ 1 , what does not change th e po- 
larization tensor (|9]) l|Misner. Thorne fc Wheelerll 19731 ). The 
gauge condition makes A — (fi — and establishes the pro- 
portionality of wa ve electric field E and A in the locally flat 
co- m oving frame (|Anile fc Breuerl Il974l ; lLandau fc Lifshitj 
Il975f ). Thus, in that frame the spatial components of a^ v 
coincide exactly with the spatial 3x3 response tensor 014k 
m (ju)i — «ifc(Aj)fc- We will derive the tensor a.ih below. It 
is usually incorporated within the dielectric tensor 



s. ^ 7rc 

£ik — Oik H TTOLik- 



(13) 



The wave equation for transverse waves in terms of Eik is 



(n T 8ik — Eik) 



£1 
E 2 



= 0, 



(14) 



where the indices i,k — 1,2, so that only the transverse 
2x2 pa rt of Eik in (e 1 , e 2 , fc) basis is t aken, and n 2 = 
/uj 2 (|Swansonll2003l:IShcherbakovll2008h. The corresp on- 



dent transport equation is ( Melrose fc McPhedranlll99lh 



dEj 
ds 



2n r c 



-AsikEk 



(15) 



where Agjj. = Eik — Sik- We take the equation (|15p , its con- 
jugate, multiply correspondingly by E£ and Ek, add and 
obtain 



d{EjE' k ) 

ds 



-(au {EtEZ) - a* kl {EiEt}) 



(16) 



for n r = 1 with the observed frequency v — Ld/(2n) neglect- 
ing emission. Again, i,k,l = 1,2. Note that a.\ 2 = — Q21 
according to Onsager principle (jLandau fc Lifshitdl 19801 ) (p. 
273). By solving the definitions of the Stokes parameters (eq. 
((6]l) for (EiEl.) and substituting the result into the equation 
(|16[1 we get the transport equation for the Stokes parameters 
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(17) 



sv 



(Vv 







"PQ «/ 



with the polarization S vector (|10p by adding the emission, 
where 

a/ = Im(a22 + a\\)jv, 

qq = Im(an — Oi22)/f, 

a v = 2Re(ai 2 )/v, (18) 

pv = 2Im(ai 2 )/V, 

Pq = Re(«22 — au)/u. 



3 DERIVATION OF RESPONSE TENSOR 

The response te n sor is d e rived for thermal pla sma in 
iTrubnikovl ||l958h : iMelrosd j 19971 ); ISwansonl J2003h . How- 
ever, there is a need to recast the derivation, since we want 
to consider both signs of charge, extend the results to non- 
thermal isotropic particle distributions, and seek for exten- 
sions to non-isotropic distributions. The importance of no- 
tation consistency cannot be overemphasized, also the ori- 
entation of our coordinate axes is different from some of the 
above sources. On the way we discover the practical signif- 
icance of the response tensor: it offers expressions for fast 
evaluation of plasma absorptivities and rotativitics. 



3.1 General isotropic particle distribution 

Throughout this subsection and next subsection we employ 
vectors and tensors in (e 1 , e 2 , e 3 ) basis, but skip tildes 
to avoid clutter. Only tildes over the response tensors are 
drawn. For simplicity of notation we define the dimension- 
less momentum 



P = Vl 2 - 1. (19) 

where 7 = En/(mc 2 ) ^ 1 is the 7-factor for particles of 
energy En and mass m. Then the Boltzmann equation on 
the distribution function f(x,p) is 

% + v ■ V/ + r,e (e + - x (B + B)) • ^ = 0. (20) 
at \ c I mc 



Here the velocity vector is 
v = pc/7, 



(21) 



Bo is the background magnetic field, E and B are the wave 
electric and magnetic fields, n is the sign of the charge and 
e > 0. Let us assume a general isotropic particle distri- 
bution fo(p) instead of thermal. The wave with the phase 
exp(i(fe ■ r — tot)) causes perturbation in a form 



/1 = exp(i(fc • r - ut))fo$(p), 



(22) 



which implements the general function &(p) of momentum 
p. Note that the perturbation (1221) is small and further anal- 
ysis is valid only when \B\ <S Bo|, wh ich corresponds t o 
low radia t ion pr essure medium. Following ITrubnikovl (|1958T) ; 
ISwansonl (|2003f ) we introduce the cylindrical coordinates 
with axis along Bo, so that 



with angle <j> in (xy) plane. Then 

k ■ p — k z Pz + k±p± sin <j>. 

The Boltzmann equation (|20[) results in 



(24) 



„ i{k z p z + k±p± sin </>)c , 
— tLUfl + /1 

7 

+ ^.M + !!! (pxBo) .M = o (25) 

mc 7 mc 

after dropping the second-order terms upon substitution of 
/1 from equation (|22|) . For general isotropic fo(p) the rela- 
tion (p x (Bo + B)) ■ Vp/o = holds, and it does not hold 
for non-isotropic distributions. Only the cylindrical <f> com- 
ponent is non-zero in a triple product (p x Bo) ■ V p /i . After 
some transformations we obtain an equation on $(p) 



zw$ + 



+ 



t(k z p z + k±p± sin (f>)c 



7 



d In /o ne 



(p-E u 



<I> 



7)U)B d§ 



= 



(26) 



dj •ymc 7 

for general isotropic particle distribution, where E u is an 
amplitude such that E — E u] exp(t(k ■ r — tot)). Here the 
cyclotron frequency is 

eBo lob 

vb = 



LOB = . • :... 

mc 2tv 

We define the ratio /3 

v 
Alternatively, the equation ()26p reads 

i(a — fesin^)$ + d<&/d<f> — F 

with 



7 — n 7.Pz 



n±pi 
VP 



F 



d In /q p ■ E u 
(by Bo 



(27) 



(28) 



(29) 



(30) 



Here n z = cos 9 and n± — sin 9 assuming jfcjc — to. The 
solution is 

$(</>) = exp(— i(a(/>+b cos (f>)) / exp(i(aip+bcos i/j))F(ip)dip, (31) 

where the lower boundary </>q is chosen at t = —00. The 
negative charge moves in the positive <j> direction and the 
positive charge in the negative cj> direction, therefore 4>o = 
—r)oo. The solution of a homogeneous equation vanishes over 
the finite time. Knowing the particle distribution and the 
definition of a current density j — ne J f\vd 3 p we calculate 
the current density amplitude 



J w = »7e / fo^(p)—d 3 p. 



(32) 



Then we relate it to the electric field E u and the vector 
potential A u wave amplitudes as 



(ju. 



Oik(E^ 



OLik{A u 



(33) 



where oiik — tuicrik/c. Let us calculate aifc. Substituting the 
solution (|31[1 into the definition of the current density equa- 
tion (I32p and changing the integration variable as ip — cf> — £ 
we get 



Px — p± cos ( 



Pv = P± sin < 



(23) 



KJw)* = t jap 



exp(— ia^ + ib(cos(<f> — f ) — cos <j))) 
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iPk{A u )k) d 



dfo 



<£. 



One can (JTrubnikovl 19581 : lLandau fc Lifshitslll980l : ISwansoij 
120031 ) introduce the differentiation with respect to vectors s 
and s' to eiiminate the momenta p in the integral expres- 
sions, then due to uniform convergence of the integrals in d£ 
and dp move the derivatives outside the integrals. We also 
do the transformation £ — > —vfi£, to finally get in (e , e 2 , 



(34) 



Qjfe = 



with 
h x ■ 



mcj 
<±me 2 d 2 



mc dsids'. 



d 2 exp(«^7 -ihp) df 



d 1 -/ 



dsids 



dfo 



jdj 



d£ (35) 



d") 



sin(ftp) 

ex P«7) 7 <% 



(1 - cos(/?£» + i(s x + cos{f3£)s' x + rjsm{j3£)s'y) 



l>v 



n± 

sm(/3£)n± 



P 



+ A s v ~ VsMPQs'x + cos(/3£)s' y ), 



h z = £n z +i(s z + s' z ), 
h = y/h£ + h v + h? z . 



(36) 



3.2 Thermal particle distribution 

Let us now consider the special case of an isotropic thermal 
distribution of particles 



/o = 



n e exp(-7/6» e ) 



4tt 9.K 2 {97 1 ) 
normalized as 



(37) 



f Q 4-Kp dp — n e 



(38) 



Here and below K n (x) is a Bessel function of the second 
type of order n with argument x. The normalized particle 
temperature is 

Q^ k A- (39) 



d£cTp (40) 



Then the response tensor (|35l) is 

d 2 exp(-A'-y -ih- p) 



Q-ik = 



„ d Sl ds' k 



^ 1 9 2 k 2 {97 1 ) 



with A' = 1/G S — it; and the rest of quantities defined by 
equation (1361) . The integral over d 3 p can be taken analyti- 
cally (|Trubnikovlll958h to give 



dik 



ie 2 n e 



d 2 



mc9 2 K 2 (9- 1 ) dsidsi 



y/A 12 + h^~ 



dg. (41) 



Performing the differentiation in Mathematica 7 to avoid 
errors, one gets 3x3 response tensor 



2*£ 



Qjfe = 

with 

71 = 

Here 

(cos /3£ t) sin /3£ 

-T^sin/JC cos/?£ 

1 



delist -*&*><*> 



7e 2 



^ 3 



e 2 sin 2 6> + 



2 sin 2 6> 



(l-COB/30- 



(42) 



(43) 



(44) 



Tik — 



-(l-cos/?£) 2 sin 2 6> 
-77(l-cos/3£)sin/3£sin 2 6> 
-VK(1- cos ££) cos 9 sing 



J 



J7(l-cos/3£)sin/3£sin 2 6l 
sin 2 /3£ sin 2 
/3£ sin /3£ sin cos 



r 



77/3^(1 — cos/3£) cos 9 sin f 
/S£ sin /3£ sin cos 9 
P 2 ^ 2 cos 2 6> 



(45) 



The expressions (I42M5[) are hiding inside two almost 
transverse and one almost longitudinal damped eigenwaves. 



3.3 Rotation of thermal response tensor 

Let us apply the transformation @ to tensors T} k and T 2 k 
and take the transverse 2x2 part to obtain correspondingly 



T 1 



cos j3£ r/ sin /?£ cos 9 

— ?7 sin /3£ cos (9 cos /?£ cos + sin 9 



and T 2 = i^i2j with 

R, = sin6>(^(l - cos/3£),cos#(sin/3£ _ ££)), 
Rj = sin0(-?j(l - cos#O,cos0(sin/3£ - /9£)) 

for 



mc9 2 K 2 {97 1 ) 



-R 2 



K 3 (TZ) 

pm 3 



(46) 



(47) 



.(48) 



The integration over £ converges very slowly, if performed 
along the real axis. The way to accelerate the convergence is 



to perform the integration in a complex plane at a positive 
angle to the real axis. The wave frequency v in the above 
calculations has a small positive imaginary part Im(z/) > 
to account for the energy pumped into particles from pass- 
ing waves. Then Im(/3) < and the expression (|43|) has zeros 
only in the lower plane Im(£) < of £ . Thus, deforming the 
integration contour to the upper plane of £ does not change 
the response tensor (|48|) . Note that all absorptivities ai, olq, 
and olv and rotativities pQ and pv are positive for electrons 
for 9 £ (0, 7r/2) under the definitions (|18[) , what gives an easy 
way to check the implementation of radiative transfer algo- 
rithm. The evaluation of these coefficients will b e repo rted 
in the subsequent paper (IHuang fc Shcherbakovl (120100 . in 
prep.). We will also evaluate the validity of a transverse ap- 
proximation for waves. 

Following IHuang et al.l (|2009T ) we define the parallel 
transported vectors a and b in addition to the right triad e 1 , 
e 2 , k, so that in the co-moving locally-flat reference frame 



(a,b) 



cos x sin x 
— sin x c °s X 



(49) 
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Then the transformation with — 2x angle 

10 

cos(2 X ) -sin(2x) 



a Q (A = 0) 



fflO) 



a a a o 



1, k"V a a a = 0, 



(59) 



R(x) = 



sin(2x) cos(2 X ) 
1 



(50) 



serves to get the vector of emissivities e and the matrix of 
rotativities/absorptivities K in (a, b, k) basis as 



R(x) 



K = R(X) 



£/ 

£Q 



olq ai 



(51) 



QV 

pv 



— pv a/ pq 
ctv — pq ai 

Define the perpendicular magnetic field 

Bo-L = B Q -k(k-B )/k 2 . 



R(-x)- 



(52) 



The trigonometric factors are related to the magnetic field 

as 

sinx = (a- B ±)/B ±, cosx = -(b ■ B ±)/B ±. (53) 
The radiative transfer equation is then 
dS/ds = e KS (54) 

for the polarization vector S defined in (a, b, k) basis. 

4 EXTENSION TO GENERAL RELATIVITY 

Let us consider two reference frames: locally-flat co-moving 
reference frame with 4- velocity u a — (1, 0, 0, 0) and flat met- 
ric and the lab frame with Kerr metric and the fluid moving 
at u a . We denote by hats () the quantities in the co-moving 
frame. Consider the radiative transfer equation ()54|) in the 
co-moving frame. The set of Stokes parameters S can be 
generalized to the corresponding set of photon occupation 
numbers 



TV = S/v 3 



(55) 



which are invariant under the orthogonal c oor dinate 
transfor mations (IMisner. Thorne fc Wheeler! 1 19731 ; 

Breued 1 1974 ; lEllid 



lAnile 



2009). 



Photons propagate 
along null-geodesies with the affine parameter A, so that 
the wave four- vector is 
dx a 



k - 



<i\ 



and 
dN 
~d\ 



0. 



(56) 



(57) 



(58) 



Here fco is a constant photon energy, which relates to the 
observed frequency as 

k c 

Under such normalization of A approximately ds « dX far 
from the BH. One calculates the null geodesic starting from 
the observer's plane. The perpendicular unit vector a a has 
a special orientation on that plane, and is transported along 
the geodesic according to 



where V CT is the covariant derivative. The unit vector b a is 
transported the same way. 

Just as in a flat space case the charged particles lead to 
the increase of occupation numbers TV due to emission, to 
decrease of TV due to absorption, and to exchange of TV com- 
ponents due to Faraday rotation and Faraday conversion. 
These are all the processes occuring in linear regime. The 
invariant number of photons emitted per unit solid angle per 
unit frequency per unit volume per unit time is proportional 
to th e invariant e{v)/v 2 = i{v)/v 2 (jMihalas fc Mihalasl 
1 19841 ) as 



riTV e{v) 

~d\ cc ~^ 

where 

v — —k^Uu 



(60) 



(61) 



is the photon frequency in the lab frame for ( — ,+,+,+) 
signature of metric. By v the photon frequency in the co- 
moving frame is denoted. The invariant change of photon 
states due to absorption and propagation effects is propor- 
tional to the co- moving frame matrix K (see eq. (l54[0 taken 
within unit frequency unit solid angle unit volume unit time. 
Thus the pr oportionality to the invar iant i/K(v) = i>K(i>) is 
established JMihalas fc Mihalasl l"l 9841 ') as 



riTV 
"dA" 



oc -^K(z^)TV, 



(62) 



where the whole form of the absorption/state change matrix 
is preserved. The full GR radiative transfer equation 



dN 
'~d\ 






(63) 



is obtained. The equation (I63p is similar t o the GR polar- 
ized transfer equation in lHuang et al.l l|2009l) . However, their 
usage of primed and unprimed quantities is potentially con- 
fusing. It is their primed quantity S' , which should be gen- 
eralized to GR as TVs = S'/u 3 . 



4.1 Transformation to locally- flat co- moving 
frame 

The angle x i n the expression (I53f) . 6 in the response ten- 
sor, and similar quantities need to either be evaluated in the 
locally-flat reference frame, where fluid is at rest, or prop- 
erly calculated in GR. We choose the first path as a trans- 
parent one with the following recipe. First, one traces the 
null geodesic from the observer's plane to the BH horizon 
or the sphere far from the BH and finds the vectors k a and 
a a (see eqs. (|56l59[) '). At each point on the ray one knows the 
vectors k a , a a , fluid four- velo city u a , and the four- vector o f 
magnetic field Bq defined in iMcKinnev fc Gammiel (2004) 
in the lab frame. The next step is to transform all vectors to 
the co-moving frame. Let us construct an orthonormal basis 
in Kerr metric in lab frame 



6(t) = (U ,U ,U ,U V ), 

e" r ) oc (u u r , — (u ut + u u,/,), 0, u u r ), 

e/$) oc (w ue, u ug,u ug + 1, u Ug), 

e (4>) K {-114,/ut, 0,0,1), 



(64) 
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where lower-index velocity is u a 
lower index Kerr metric 



9afi 



( -1+' 




2ar sin 



V- 2 





p/A 







a pw and (? a /3 is the 
\ 



lar sin 2 6, 







iin 2 e„ 



(65) 



/ 



in (t,r,6 a ,4>) spherical polar coordinates with polar angle 
8 a , radius r, spin a, p = r 2 + a 2 cos 2 6 a , A = r 2 — 2r + a 2 , 



E = 
to 



(r 2 + a 2 ) 2 



a A sin 9 a . Then make a transformation 



ef t) = (1,0,0,0), 
Ig.) = (0,1, 0,0), 

ef e) = (0,0,1,0), 
e^ = (0,0,0,l) 

via 



(66) 



)(t,r,9,<j>)p 



(-e( t) ,e 



(t),e (r ),e (e) ,e (0) 



U))9a0- 



(67) 



The transformation of a four-vector A^ to the co-moving 
frame is then 



The metric in the new frame 



(68) 



(69) 



coincides with Minkowski metric rjik = diag(— 1, 1, 1, 1). The 
velocity four-vector u 13 transforms to um = (1,0,0,0) T . 
Thus, the basis change (168[) with matrix (|67|) and vec- 
tors (|64|l constitutes the transformation to the locally-flat 
co-moving reference frame. This procedure is the alter- 
native of the numerical Gramm-Schmidt orthonormaliza- 
tion applied in iMoscibrodzka et al.l <I2009| N ) . The ba s is vec - 
tors (164)) are presented in Krolik. Hawlev fc Hirosd (J2005J); 
iBeckwith. Hawlev fc Krolikl ( 20081). our expressions being a 
simpli fied version of vectors in IBeckwith. Hawlev fc Krolikl 
(2008). Note, that despite the vectors (|64[) do not explicitly 
depend on the metric elements, the expressions rely on the 
properties of Kerr metric and are not valid for general g . 
Upon transforming u a , k a , a a , Bq to the co-moving 
frame we easily find the wave frequency v = — k(o)> then 
k — fc(i,2,3)- The perpendicular vector a needs to be offset 
by k as 



,- a (o) 

a = 0(1,2,3) - k- 



k. 



(70) 



(") 



and then normalized to construct a spatial unit vector. The 
offset is due to the enforcement of Lorenz gauge <|12[) . It 
conveniently makes a ■ k = a ■ k = 0. The vector b is found 
by a simple vector product 



a x k 



(71) 



and then normalized. The spatial part (e( r ), e($) : &(<!>)) of 
basis (|66[) relates to basis (e , e 2 , e 3 ) via the orthonor- 
mal transformation preserving angles. The magnetic field 
Bq gets transformed to a three-vector Bo and all the an- 
gles are found in correspondence to Fig. \T\ with the help of 
equations (J52I53I) applied to hatted vectors. For example, 
cos# = (k ■ Bo)/(kBo). Then the whole matrix K is found. 



5 APPLICATION TO COMPACT OBJECTS 

The described GR polarized radiative transfer finds its ap- 
plication in accretion onto low luminosity active galactic nu- 
clei (LLAGNs), in particular Sgr A*. The application of GR 
is necessary to infer the BH spin, which provides impor- 
tant information on the past evolution of the BH and the 
host galaxy itself. For example, the accretion efficiency in 
the AGN pha se depends strongly on the value of BH spin 
(Shapiro 2005). The value of spin and spi n orientation con- 
strain ts the accretion and merger history (|Rees fc Volonteril 
l2007l V 

The detailed a pplicatio n to Sgr A * is reported in 
IShcherbakov. Penna fc McKinnevI (|2010f ) . Let us describe 
on the example of that paper, how one connects to obser- 
vations. First, one constructs a set of dynamical models, 
preferably based on 3D GRMHD simulations. Then the GR 
polarized radiative transfer is performed for those models as 
described in the present work. The simulated spectral energy 
distributions (SEDs), linear (LP) and circular (CP) polariza- 
tion fractions as functions of frequency foo are fitted to the 
observed quantities, representative for the quasi-quiescent 
state of accretion. The x 2 analysis is performed based on 
the inferred error or variability of the observed quantities. 
Then the best fits in the parameter space can be found. The 
probability density can be integrated over the full parame- 
ter space to obtain the most likely values and the confidence 
intervals of BH spin, inclination angle, position angle, and 
model parameters. 

Let us now describe the effects, which lead to certain 
observed cyclo-synchrotron spectra, LP and CP fractions, 
and the electric vector position angle (EVPA) as functions 
of frequency. The effects are plenty, what proves it hard to 
disentangle and provide simple explanations of observations. 
It is in general challenging to achieve the realistic level of 
details in collisionle ss plasma mo d elling . The next step in 
IShcherbakov. Penna fc McKinnevI (]2010r ) might not be the 
most self-consistent. First, the radiation from LLAGNs ap- 
pears to be variable in time. The simultaneous short obser- 
vations can provide only the single snapshot of a system, not 
necessarily representative of a long term behavior. Thus, it 
is necessary to obtain the statistically significant sample of 
variability of both observed and simulated spectra to reli- 
ably estimate the average or typical flow parameters and BH 
spin. As recent research suggests (Dodds-Edcn ct al. 2010), 
the modelling of a single flare can successfully be done even 
without invoking GR. 

As in the case of Sgr A*, the cyclo-synchrotron specific 
flux F v vs v can have a peak. The peak frequency v* and flux 
F£ do not necessarily correspond to the thermal cut-off of 
emission. Even a small percentage of the non-thermal elec- 
trons can radiate significantly more than the bulk of ther- 
mal electrons. For the efficient particle acceleration most of 
emission may come from the energetic electrons with cooling 
time t coo ] about the time of inflow £j n from several BH grav- 
itational radii r g = GM/c 2 . The gravitational redshift and 
Doppler shift due to relativistic motion can strongly modify 
the peak v* and F* . 

The LP fraction can provide constraints on flow den- 
sity near the emitting region owing to beam depolarization 
effect. The LP fraction is the highest at high frequencies, 
where only a small region of the flow shines and beam de- 
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polarization is weak. As all regions of the flow radiate at 
lower frequencies, differential Faraday rotation and emission 
EVPA vary, the resultant LP fraction is subject to cance- 
lations and quickly ceases with v. However, cancelations at 
high v may readily happen between two regions with similar 
fluxes and perpendicular EVPAs, those regions would have 
perpendicular magnetic fields. 

The same change of EVPA with frequency can mimic 
the Faraday rotation. The finite rotation measure (RM) 



RM 



EVPA! - EVPA 2 
X 1 - A 2 



(72) 



does not necessarily happen due to Faraday rotation ~ 
J nB ■ dl. Here A = cjv is the wavelength. In fact, the 
meaningful application of formula (|72[) is limited to a toy 
case of cold plasma far from the emitting region with the 
homogeneous magnetic field. In reality, besides the change 
of emission EVPA with v, the Faraday rotati on coefficient 
Pv(v ) (see eq |18|) is a function of frequency IjShcherbakovl 
I2OO81 ). The differential rotation measure 

dRM = dEVPA/d(\ 2 ) (73) 

is the measured quantity (JMarrone et al.l 120071 ). It should 
be used in constraining the models. Also, significant Fara- 
day rotation can happen in the emitting region, what in- 
troduces the effect of differential optical depth. Thus, one 
can only fit the observed EVPA(^) and use it along with 
other observables to constrain the system free parameters. 
As py(v) is a steeply declining function of temperature 6 e 
(Shcherbakov 2008), the relativistic charges contribute little 
to this quantity. 

Substantial l evels of circul a r pola rization were recently 
found in Sgr A* (|Munoz et al.l (|2010l '). in prep.). There are 
several effects producing finite CP. First recognized was the 
emissivity ev in V mode. According to Melrose (1971) it 
only is a factor of 7 weaker than the total emissivity ej. 
It produces the largest V along the magnetic field. The 
Faraday conversion, transformation between the linear po- 
larization and the circular, operates perpendicular to the 
magnetic field. The Faraday conversion coefficient pQ has 
a peculiar dependenc e on temperature o f thermal plasma 
or particles' 7-factor l|Shcherbakovll200a ): pq = for cold 
plasma, pq is exponentially inhibited for very hot plasma 
and reaches the maximum for transrelativistic plasma with 
6 e ~ 7 ~ several • mc 2 /fe. The exponential inhibition is an 
effect of finite ratio v /vb with peak pQ only around 8 e ~ 10 
for vjvs ~ lO 3 . Thus, for hot bulk part of particle dis- 
tribution with 7 ~ several • mc 2 jku the non-thermal elec- 
trons do not contribute significa ntly to Faraday conversion. 
Note that this result supersede s Ballantvne. Ozel fe Psaltisl 
l|2007l ). who following iMelrosd (|l997h neglected the impor- 
tance of finite ratio v/vb- Similarly to linear polarization, 
the beam depolarization can lower the net value of circular 
polarization at low frequencies due to differential Faraday 
conversion. 

In sum, there are often several explanations for the same 
observable quantity. One should not settle for a simplified 
model trying to reproduce the observations. Instead, the rig- 
orous ray tracing and aposteriori explanations of a fit to the 
observables would be the preferred reliable way. 

The proposed method has its limitation. The equation 
(|63|) is valid for optically thick medium, but it fails to de- 



scribe the behavior of a set of photons for Compton-thick 
medium. The encounters of photons with energetic electrons 
lead to significant changes in photon trajectory, whereas the 
previous discussion considered independent photons prop- 
agating along geodesies. Luckily, the synchrotron absorp- 
tion cross-section in sub- mm is much larger than the Comp- 
ton scattering cross-section, thus the optically thick medium 
near Sgr A* is Compton-thin and no modifications are 
needed for Sgr A*. However, a careful consid eration of 
Compton scattering (jRvbicki fe Lightmanlll979l ) is needed 
to describe the sub- mm spectrum of Compton-thick sources. 



6 DISCUSSION & CONCLUSIONS 

In our endeavor to provide the complete and self-consistent 
description of GR polarized radiative transfer we conduct 
the full derivation starting from definitions and basic equa- 
tions. The goal is to make the easy, transparent, and error- 
free derivations, thus Mathematica 7 was used underway, 
the expressions were cross-checked. The absorptivities for 
thermal plasma were checked numerically against known 
synchrot ron emis si vities and cyclo-synchrotron a pproxima- 
tions in ISazonovl <|l969h : iLeung. Gammie fc Nobll (|2009h . 
We stepped away from the standard textbooks and assumed 
"the opposite" observers' definition of circular polarization 
V, carrying the definition through all other calculations. We 
chose the coordinate system with coplanar fc, e 2 , Bo and 
derived the plasma response tensor &ik in (e 1 , e 2 , e 3 ) basis, 
also projecting it onto th e transverse coordin ates e 1 , e 2 . Re- 
peating for completeness iHuang et al.l (|2009| ) , we tie the po- 
larized radiative transfer equation in the latter coordinates 
(|17|) with the transfer in a, b coordinates with the help of 
matrix (|50[) . The generalization of radiative transfer to GR 
is performed in the easiest way owing to the invariance of oc- 
cupation numbers in different photon states and the invari- 
ance of the transformation between the states. Lorenz gauge 
H12|) helps to establish the correspondence between 4x4 and 
spatial 3x3 quantities and to correctly find the transverse 
vectors a, b in the co-moving locally-flat reference frame. 
The transformation from the lab frame with Kerr metric to 
that frame is explicitly given. The intricacies of application 
of GR polarized radiative transfer to LLAGNs are discussed. 
As the transfer incorporates many physical effects, a priori 
guessing of the most important effects is discouraged in fa- 
vor of full calculation. The provided interface of dynamical 
models and observations is waiting for its applications. 

The treatment of particles distributions is still lim- 
ited. In the current state the calculations are optimized for 
isotropic in pitch-angle distributions and become especially 
simple for thermal particle distribution. The integration over 
the pitch-angle in formula (|35|) is in general impossible to 
perform for non-isotropic distributions. In this case, the inte- 
gral over £ should be done first analytically. This calculation 
is left for future work. 
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